Energy and entropy of metastable states in glassy systems 
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We investigate the multi-valley energy landscape of a 3-D 
on-lattice network model for covalent glasses, numerically de- 
termining the shape of the valleys, the local density of states, 
the density of minima and the local connectivity. We present 
some of these quantities in a graphical birds-eye view of the 
landscape, and discuss their implications for the relaxation 
dynamics and cooling behavior of glasses. The strong sim- 
ilarities between the landscape of this model and those of 
other complex systems point to the possibility of a common 
low-temperature dynamical description. 
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Introduction. A number of problems ranging from 
solid state physics to combinatorial optimization [Q, 
through chemical physics 1^-^, and biology ^,0] can 
be formulated in terms of Marcov processes with states 
in the systems' energy landscapes. The (effective) con- 
nectivity of the landscape specifies the allowed hopping 
transitions while the energies determine their magnitude. 
Specifically for the multi- valley landscapes characterizing 
complex systems, local thermal equilibrium establishes 
itself within 'ergodic components' and an investiga- 
tion of the local geometry close to energy minima yields 
significant dynamical information. 

For discrete systems, all configurations within subsets 
of the state space (pockets) centered around local energy 
minima can be counted |9|. In addition, combinatorial 
techniques have been used to find all low-lying minima 
of small systems , and different approaches have been 
devised to coarse grain the configuration space into a 
space of more manageable size [^p|Jll|-[l^ . 

For clusters and glassy systems with continuous config- 
uration spaces the main emphasis has been on stochastic 
methods for finding energy minima and the lowest sad- 
dles which connect them Q. The vibrational modes of 
the systems can then be described within a harmonic ap- 
proximation of the potential close to the minima. 

Low temperature MC/MD-simulations have also been 
utilized to probe the vibrational and relaxational be- 
haviour on short time scales in the neighborhood of 
" frozen-in- networks" |14-1^]. At higher temperatures, 
these methods yield insights into the dynamics of the 
glass transition and aging phenomena [f? 18|. In addi- 



tion, the large scale distribution of local minima 20 
can be sampled by quenching out of the melt. 

While knowing the minima and their connecting sad- 
dles is useful in predicting e.g. whether a given system 
will be glass forming or not quantitative knowl- 

edge of the density of states of model glasses is still rather 
scanty. This information can be useful in clarifying the 
connection between structural glasses and systems with 
quenched randomness as e.g. spin-glasses and help to 
provide a unified picture of the energy landscapes of com- 
plex systems. 

With this purpose in mind we present in this Let- 
ter large-scale exhaustive enumeration of microscopic 
states for on-lattice 3-D network models of glassy systems 
pl|-p3|. By placing the system on a lattice, a procedure 
successfully employed for polymers |^,^, we can ana- 
lyze in detail regions of configuration space containing 
many thousands of local minima, each minimum repre- 
senting a topologically different network. 

Our analysis shows that the local trap geometry can 
be characterized rather simply. The discussion deals with 
the connection between the geometry and dynamical fea- 
tures of glasses, including thermal metastability, domain 
formation and cooling behavior. 

Model. The main structural feature of covalently 
bonded glasses as e.g. amorphous Si02 close to the tran- 
sition temperature is their topology. Accordingly, such 
systems are often modeled as random networks of build- 
ing blocks where links represent the covalent bonds. 

Since we are interested in configurational changes 
which take place on time and length scales considerably 
larger than those appropriate for the vibrational degrees 
of freedom of the system, we can put our network on a 
grid with a lattice spacing a of the order of the bond 
lenght. The harmonic vibrations cannot be described at 
this level of resolution. Furhermore, the energy minima 
of the discrete model and the corresponding positions of 
the atoms will be somewhat changed in comparison with 
the off-lattice network. However, overall geomerical fea- 
tures of the energy landscape as e.g. the distribution of 
accessible minima and barriers are not likely to be seri- 
ously affected. 

We use cubic lattices with lattice constant a = 1, pe- 
riodic boundary conditions and a repeated cell of linear 
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size S = 10 and volume V = S^. The number of atoms 
is Na and the density is p = Na/V. 

There is an arbitrariness in the choice of the elementary 
moves which define the dynamics on the energy landscape 
and hence its connectivity. We were guided by the simple 
physical consideration that a single move should involve 
a change of coordinates which is small. Each elementary 
move is therefore taken as a shift of the position of one 
atom to one of the six adjacent grid points. With this 
definition each state space configuration has 6Na neigh- 
bors in the landscape. 

Our interaction potential is constructed such that a 
well defined crystalline ground state is present at a cer- 
tain density. At the same time it allows for many other 
stable configurations with no crystalline structure. The 
potential vanishes except for nearest neighbors, which 
are atoms closer than r = 3.5. The nearest neighbor in- 
teraction consists of a two-body and a three-body term. 
The radial dependence of the former has a hard-core for 
r < 1.6, and equals V2(r) = (r - 1.97)(r - 3.5)2(r - 7), for 
1.6 < r < 3.5. This potential has a minimum at r « 2.2. 

The three-body term ¥3(9) containing the angular de- 
pendence is infinite for angles below 40°. In the range 
[40°, 80°] it is strongly repulsive and given by ¥3(6) = 
-2.9 X 10-^(6* - 80) (61 - 180)26*. In the range [80°, 180°] 
the potential has a shallow broad minimum close to 120° 
and vanishes at the high end of the interval, as expressed 
by the form: ^3(6!) = -5.8 x 10-^(6' - 80)(6' - 180)26*. 
The optimal local coordination in the dense phase is the 
result of the interplay of V2 (r) and V3 (9) , the former fa- 
voring high coordination (up to six neighbors in an octa- 
hedral arrangement), and the latter one favoring a trig- 
onal planar environment. As a consequence, four-fold 
(tetrahedral) coordination is favored energetically, fol- 
lowed closely by five-fold coordination. Configurations 
connected by translations and/or rotations of the sys- 
tem and hence physically equivalent are identified and 
counted once. 

Results. The parallel version of the lid-algorithm uti- 
lized for our exhaustive searches is detailed in Ref. [^ . 
The basic idea of how the method works is best explained 
by referring to the data shown in Fig. 1. For concrete- 
ness we shall focus on the cone-like structures pertaining 
to the 50 atom system. The algorithm starts at the ref- 
erence configuration represented by the apex of the high- 
est lying cone and visits all configurations which can be 
reached by a succession of elementary moves never ex- 
ceeding an energy barrier L. The latter appears as the 
vertical distance from the apex along the cones' symme- 
try axis. The volume V{L) of a pocket of depth L is then 
the number of configurations reachable below energy L. 
In the figure, the width of the cone at distance L from the 
apex is proportional to the logarithm of V{L). As we see, 
the accessible volume is almost exponentially increasing 
with the barrier. At some value Lmax of the barrier, 
which in the case of the highest pocket is ~ 0.5eV, the 



algorithm finds a state of energy lower than the lowest 
energy so far. This state is adopted as a new reference 
state, and the counting starts afresh. By iterating this 
procedure the algorithm generates a succession of cones 
located at progressively lower energies. Each of these 
is the graphical rendering of a state-space pocket. We 
emphasize that, while the counting within each pocket 
is exhaustive, the path joining one pocket to the next 
(symbolized by the lines connecting the cones) is just one 
amongst the many possible. The largest pockets contain 
w 10^ — 10^ states. A second quantity of importance is 
the local density of states ^{E; L), describing the energy 
distribution of the states in a pocket of depth L. For 
some representative instances InV is plotted versus E, 
where < < L, in the same way as V{L) was plot- 
ted versus L. This results in the conical shapes enclosed 
in the dashed boxes of Fig. 1, which show that ^{E; L) 
is also close to an exponential function of its argument. 
A more detailed analysis shows that the curves have a 
slightly concave appearance in a semilog plot. 

Figure 1 may convey the false impression that the 
pockets are devoid of internal structure. In reality, they 
all contain a multitude of local minima, each surrounded 
by its own basin. In fact, the number of accessible min- 
ima and their density also grow approximately exponen- 
tially as a function of L and E (for fixed L) , respectively. 
Part of this rich structure is made visible in Fig. 2, where 
we consider one of the previously depicted pockets (deep- 
est minimum for Na = 45) in more detail. The local den- 
sity of states of a pocket (top curve) and six of its sub- 
pockets are displayed. We remark that each sub-pocket 
contains many local minima, and all the local densities 
of states have a similar, close to exponential appearance. 
We have also investigated the distribution in energy of 
the local minima in each pocket finding shapes similar 
to that of the local density of states, except for stronger 
flattening close to the lid energy. 

In a harmonic approximation of the vibrational modes 
within a local minimum, one would expect a power-law 
form of the density of states. The present local densities 
of states, however, pertain to a rather different situation, 
i.e. the much larger length and time scales of excitations 
associated with changes in bonding topology. To check 
whether a power-law could nonetheless describe the data, 
we shift the abscissa with the energy of the lowest min- 
imum in each valley and plot the data in a log- log plot. 
The inset of Fig. 2 shows two of the curves of the main 
panel, plotted in the fashion just described. 

Finally, we considered the radial and angular distri- 
bution functions G(r) and P{9) for networks of different 
densities. In Fig. 3, we depict P{9) averaged over the 
lowest minima for each density. The fact that the strain 
in the network increases with density is revealed by the 
increasing predominance of 'bad' angles. 

The acquisition process for all of our data took of the 
order of one year of calculation, running in parallel on 8 
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processors on a SGI Onyx computer. 

Discussion. The local geometrical structure of each 
pocket can be described in terms of available volume, lo- 
cal density of states, and local density of minima. Even 
though these quantities vary somewhat from pocket to 
pocket, they nevertheless share important features. Ap- 
proximately, the local density of states within a pocket of 
depth L is V{E,L) k. C{L) e^p{E/Tg). The parameter 
1/Tg characterizes the average slope of I? in a semilog 
plot and the function C{L) describes the effects of side- 
valleys joining the main valley as the lid grows. Inte- 
grating V with respect to E from E — up to E ^ L 
yields V{L) oc C{L)exp{L/Tg). Since, as seen in Fig. 1, 
V{L) also grows almost exponentially as a function of 
L, we conclude that C{L) must itself be an exponen- 
tial or equal to a constant. In the latter case all the 
states added by increasing the lid would have energies 
very close to the lid itself, which in turn implies that the 
local minima are shallow or lacking altogether. In fact, 
however, deep sub-pockets frequently appear as the lid 
is increased. Since these side-pockets have fewer states 
than the main pocket, but the same overall shape with 
a lesser barrier to a deeper minimum, they are approxi- 
mately scaled-down versions of their parent. Referring to 
Fig. 1 we note that this scaling property also seems to ap- 
ply when moving one level up from the individual pockets 
to the larger scale structure of the whole landscape. This 
is because a high-lying pocket is by construction a sub- 
set of a lower-lying one and because all pockets within 
a certain energy range seem to have similar shapes inde- 
pendently of the energy of their lowest minimum. All in 
all, this is strongly suggestive of a self-similar structure 
in the energy landscape. Of course, to rigorously check 
self-similarity one should extend the analysis to several 
additional levels of nesting, a task which falls outside the 
scope of this investigation. From what we can glean so far 
concerning the largest energy scale, where we also would 
include very high-lying (possibly physically non-relevant) 
minima, there appears to be a weak trend to higher values 
of Tg for the high-lying pockets. But since these narrow 
pockets usually appear to be relatively small compared to 
the rest of the system, their influence on the dynamical 
properties should not be very strong. Quite generally, the 
information gained from these exhaustive searches allows 
the construction or/and validation of mesoscopic models 
of complex system dynamics p8|-|30| which rely on as- 
sumed geometrical properties of state space. Since our 
pockets are metastable below the trapping temperature 
Tg one can construct and use the connectivity matrix of 
small pockets (< 10^ states) to obtin direct insight in the 
low-temperature relaxation behavior ^ . 

For T > Tg, a hypothetical quasi-equilibrium distribu- 
tion within the pocket would be very skewed in favor 
of states with energy close to the lid [p pl|j3^ . Ban- 
ning unlikely kinetic mechanisms which would prevent 
the system from escaping the trap, this implies a loss of 



metastability at T « Tg, where the trap effectively disap- 
pears from the landscape. A purely exponential density 
of states gives a singularity in the specific heat when the 
temperature approaches Tg from the low side. (The high 
side is irrelevant since the system has then left the pocket 
and new unchartered parts of state space must enter the 
calculation). Such a peak in Cy is also observed exper- 
imentally at the glass transition temperature, and com- 
monly associated with a configurational entropy pl| , p2[ . 
The finite size of the pocket and corrections to the purely 
exponential behavior smoothen the transition of course. 

The value of Tg decreases with the system real space 
volume V ||2^ , |3l| . In the glass case, a simple free volume 
analysis |23[ indicates that Tg oc 1/(1^(1 — pVa)), where 
Va is the volume/atom in a dense strain-free network 
and p is the particle density. The fact that Tg vanishes 
in the limit of infinite volume does not imply that trap- 
ping does not take place in large systems. Indeed, a large 
system will not be able to equilibrate internally at low 
temperatures, and the interpretation of Tg{V) as a trap- 
ping temperature is then no longer appropriate. Instead 
one may envisage that a large system may form domains 
of strongly correlated atoms, similarly to the behavior 
observed in spin-glasses ]3^|35[. The preceding analy- 
sis applied to these domains leads to the conclusion that 
the domain size is limited by an entropic mechanism: for 
a given temperature T a domain can only be stable if 
Tg{Veff) < T. Domains can grow up to the size at which 
equality holds in the above relation. 

If the system is annealed from high temperatures a 
slower cooling rate will arguably lead to a larger size 
Vef f of locally equilibrated domains, and hence to a lower 
glass transition temperature, which agrees with the ex- 
perimental behavior of glasses under cooling |2l|,|2^ . 

The basic condition for the above scenario, i.e. the ap- 
proximate exponential growth of local densities of states 
in locally ergodic regions of the energy landscape, is ful- 
filled in a number of additional systems that show com- 
plex dynamics and glass-like behavior: the TSP spin 
glasses Q, and lattice polymers Thus, we suggest 
that a loss of metastability in spatially localized traps 
with exponential densities of states, i.e. mainly an en- 
tropic mechanism, is an important factor in the glass 
transition and in cooling rate dependent dynamic effects. 
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valley index (arb. units) 



FIG. 1. A birds-eye view of the energy landscape of the 
glass model, where pockets appear as cone-like structures. 
The width of the cones at energy L relative to the position 
of the apex is proportional to the logarithm of the locally ac- 
cessible state-space volume V{L). The insets use the same 
graphical rendering but with the cone-width at energy E pro- 
portional to the logarithm of 'D[E\ Lmax), the local density of 
states available in the selected pockets. 
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FIG. 2. The local density of states of several nested 
sub-pockets of a system with 45 atoms. The top curve (black 
squares) describes the local density of all states in the full 
pocket up to an energy of « 0.22 eV/atom. Each other curve 
describes the local density of states of a sub-pocket merging 
into the main valley approximately at the energy where the 
curve ends, within the accuraqcy given by the binning of the 
energy axis. E. g., the second curve (empty triangles) corre- 
sponds to a subpocket merging at « 0.19 eV/atom. Note that 
the growth of the local densities of states of the sub-pockets 
closely parallels the one of the main pocket. The inset shows 
the same data as the top and bottom curves of the main panel, 
plotted on log-log scales and with the abscissa values shifted, 
such that the lowest energy in the appropriate (sub)valley is 
located at the origin. 
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FIG. 3. The length of each ray is proportional to the fre- 
quency of the corresponding bond angle in the configurations 
of the deepest energy minima of the three lowest pockets for 
systems of 35, 40, 45, 50 atoms, respectively. The 30 atom sys- 
tem (not shown) is similar to the 35 atom one, exhibiting even 
more (ideal) tetrahedral surroundings. The preferred bond 
angle is close to 110°, due to the interplay of the two-body- 
and the three-body-terms in the interaction potential. Note 
the increasing predominance of energetically unfavorable an- 
gles as the atom density increases. 
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